Predicting River Resilience to Disturbance Using Population Models at a Continental Scale

Purpose Statement

We are interested in learning about streams’ responses to disturbance (i.e., storm events). Using a state-space time series approach, we are fitting a logistic population growth model to gross primary production estimates to gain information about growth parameters (e.g., maximum growth rate, critical disturbance thresholds) that describe algal population dynamics and compare results across stream sites.

Revised Hypotheses

How does the maximum growth rate (r) and carrying capacity (K) vary in streams and rivers across the country as a function of stream order, disturbance, and light availability?

H1: As stream order increases, the maximum growth rate will decrease (due to size, water depth, light attenuation and sediment load in water) and the carrying capacity will increase (due to size and space availability).

H2: The maximum growth rate will be greatest in streams with regular disturbance events that scour existing biomass. Carrying capacity will be altered if the regular disturbances are of large enough scale to alter the physical structure of the stream (e.g., could decrease due to infill with sediment or increase due to increase in surface area).

H3: As light availability increases, the maximum growth rate will increase (photoinhibition will likely not be a concern except in intermittent/ephemeral and low order streams) and the carrying capacity will also increase (up until a certain point, where the demand for light is fulfilled and no space remains).

Model Structure

A conceptual diagram of the model structure used may be found in Figure 1.

Since our last meeting, we experimented with structuring the model in a hierarchical format so as to share information across sites. However, this approach ended up being too computationally intensive, so we have reverted to modeling each site individually.

Figure 1. The state-space population model being used to estimate parameters of maximum growth rate (r) and carrying capacity (K) across individual stream sites. The observation model uses daily estimates of gross primary production (g) and daily light (L) to estimate daily biomass (B). The daily persistence term (P) is a function of daily discharge (Q) which is used to estimate the sensitivity of the persistence transition from presence to removal of biomass (s) and the critical discharge threshold at which benthic biomass is removed (c). Both the observation and persistence models then inform the process model, which is a density-dependent growth function for daily biomass (B) in a given stream. The estimated terms of P and B are input into the process model to extract values for the maximum growth rate (r) and lambda, which is equal to the maximum growth rate divided by the carrying capacity (K).

Figure 1. The state-space population model being used to estimate parameters of maximum growth rate (r) and carrying capacity (K) across individual stream sites. The observation model uses daily estimates of gross primary production (g) and daily light (L) to estimate daily biomass (B). The daily persistence term (P) is a function of daily discharge (Q) which is used to estimate the sensitivity of the persistence transition from presence to removal of biomass (s) and the critical discharge threshold at which benthic biomass is removed (c). Both the observation and persistence models then inform the process model, which is a density-dependent growth function for daily biomass (B) in a given stream. The estimated terms of P and B are input into the process model to extract values for the maximum growth rate (r) and lambda, which is equal to the maximum growth rate divided by the carrying capacity (K).

Figure 3. Persistence curves for a site that performs well with the model structure described above (South Branch Potomac River) and a site that performs poorly (Reedy Creek). The critical discharge at which biomass is disturbed (c) is denoted by the dotted vertical line. The sensitivity of the persistence transition from presence to removal of biomass (s) is equal to the slope of the green persistence curve as it transitions from complete persistence (P = 1) to removal (P = 0) of biomass.

Figure 3. Persistence curves for a site that performs well with the model structure described above (South Branch Potomac River) and a site that performs poorly (Reedy Creek). The critical discharge at which biomass is disturbed (c) is denoted by the dotted vertical line. The sensitivity of the persistence transition from presence to removal of biomass (s) is equal to the slope of the green persistence curve as it transitions from complete persistence (P = 1) to removal (P = 0) of biomass.

Data Acquisition

Data have been assembled using the datasets available in Appling et al., 2018, Blaszczak et al., 2021, Bernhardt et al., 2022, and Savoy and Harvey, 2021.

We have revised the following rules to filter the available data:

  • Filter for high-quality sites based on model diagnostics (e.g., K600 sigma Rhat < 1.05, GPP < 0 or ER > 0 for more than 15% of estimates).
  • Filter for sites that have a maximum time gap of 14 days and a minimum of 90 days per year of GPP data.
  • Any GPP estimates below zero were assigned a value between 0.05 and 0.13.

From the initial 356 river sites in the Appling dataset, our filtered dataset consists of 182 sites and 691 site-years of data (Table 1).

Table 1. Summary statistics for all 182 stream sites included in the dataset.

Table 1. Summary statistics for all 182 stream sites included in the dataset.

Revised Results

Figures of all available covariates (discharge, temperature, light) as well as the full set of our model diagnostics are available to view on a separate R Shiny application located here. Below are a few overview figures examining initial relationships among model parameter values and site covariates.

Prior to plotting, all sites whose maximum growth rate parameter had an Rhat value greater than 1.05 have been removed (n = 159 sites remaining).

Figure 4. Plots of mean maximum growth rates as (A) an overall distribution across all sites and versus (B) the coefficient of variation in discharge, (C) daily mean light availability, (D) stream order, (E) site latitude, (F) site longitude, (G) watershed size, and (H) primary land use.

Figure 4. Plots of mean maximum growth rates as (A) an overall distribution across all sites and versus (B) the coefficient of variation in discharge, (C) daily mean light availability, (D) stream order, (E) site latitude, (F) site longitude, (G) watershed size, and (H) primary land use.

Prior to plotting, all sites whose critical disturbance threshold had an Rhat value greater than 1.05 have been removed (n = 141 sites remaining).

Figure 5. Plots of the ratio of the critical discharge (Qc) to the two-year flood (Q2yr) as (A) an overall distribution across all sites and versus (B) the coefficient of variation in discharge, (C) daily mean light availability, (D) stream order, (E) site latitude, (F) site longitude, (G) watershed size, and (H) primary land use.

Figure 5. Plots of the ratio of the critical discharge (Qc) to the two-year flood (Q2yr) as (A) an overall distribution across all sites and versus (B) the coefficient of variation in discharge, (C) daily mean light availability, (D) stream order, (E) site latitude, (F) site longitude, (G) watershed size, and (H) primary land use.

Next Steps: Following this meeting, we intend to pursue two main activities:

  1. Finish a full draft of the manuscript. Please see the working Google document here.

  2. Finalize model diagnostics to determine the list of sites included in the results/discussion.

Questions for You:

  1. Do you have any suggestions regarding model diagnostic best practices?

  2. Do you feel there are particular figures that should/shouldn’t be included in the main text?

  3. Is there anything else you feel I should consider examining moving forward?

End of RMarkdown document.